%PRIME THE SIMULATION
%Define matrices to map interventions across runs
fiveyear=zeros(10000,1);
tenyear=zeros(10000,1);
fifteenyear=zeros(10000,1);

fivezero=0;
fiveone=0;
fivetwo=0;
fivethree=0;
fivefour=0;
tenzero=0;
tenone=0;
tentwo=0;
tenthree=0;
tenfour=0;
fifzero=0;
fifone=0;
fiftwo=0;
fifthree=0;
fiffour=0;
Icounts=zeros(8,3);


%Input initial country variables as row vectors, in the order they appear
%in the main model
Bahrain=[1,0,0,0,1,.396,1,0,7.47,0,0,0,0,1,log(40000),0,0,0,1];
Iran=[0,1,0,0,1,1.031,0,0,7.47,0,0,0,0,1,log(57000),1,0,0,1];
Jordan=[1,0,0,0,1,.761,0,0,7.47,0,0,0,0,1,log(500000),1,0,0,1];
Kuwait=[1,0,0,0,1,.458,1,0,7.47,0,0,0,0,1,log(40000),1,0,0,1];
Qatar=[1,0,0,0,1,.411,1,0,7.47,0,0,0,0,1,log(40000),0,0,0,1];
Saudi=[1,0,0,0,1,.877,0,0,7.47,0,0,0,0,1,log(40000),1,0,0,1];
Syria=[1,0,0,0,1,.954,0,0,7.47,0,0,0,0,1,log(1200000),1,0,0,1];
Turkey=[0,0,1,0,1,1.032,0,0,7.47,1,0,0,0,1,log(10000),1,0,0,1];

%Enter coefficients as column vectors
coeffsall=[.4664, -.9752, 1.1088, -.1959, -.0646, .2307, -.4586, .1705, .0816, -.2629, .4762, .2678, .4029, .0554, .0170, 1.4917, .2696, .0351, -4.7628]';
coeffsallb=[.513,-.909,1.105,-.182,-.075,0.236,-.477,0.156,0.082,-.257,0,.295,.396,.083,.018,1.477,0.286,0,-4.747]';
coeffsprd=[.564, -.675, .766, -.121, -.178, .249, -.447, .045, .080, -.289, .354, .303, .444, .091, .012, 1.117, .269, .113, -4.310]';
coeffsprdb=[0597,-.630,.763,-.010,-.201,.255,-.453,.035,.084,-.297,0,.321,.449,.133,.013,1.102,.330,0,-4.311]';
coeffscontig=[.483, -.458, .305, .148, -.178, .407, -.345, -.115, .105, -.252, .623, .268, .444, .078, .012, .854, .174, .245, -4.334]';
coeffscontigb=[.560,-.454,.372,.149,-.262,.368,-.521,-.103,.116,-.236,.691,.232,.452,.034,.015,0,.142,.123,-3.404]';
coeffscontigc=[.544,-.392,.340,.174,-.225,.419,-.347,.131,.110,-.269,0,.306,.452,.159,.015,.851,.311,0,-4.368]';
coeffscontigd=[.536,-.434,.360,.156,-.204,.411,-.372,.137,.104,-.260,0,.308,.413,.126,.014,.871,.192,.231,-4.368]';
coeffsall2=[.600, 0, .918, 0, 0, .206, -.420, 0, .072, -.303, .458, .236, .366, 0, 0, 1.489, .250, 0, -4.568]';
coeffsprd2=[.661, 0, .530, 0, 0, .207, -.401, 0, .077, -.302, .407, .255, .413, 0, 0, 1.06, .273, 0, -4.187]';
coeffscontig2=[.572, 0, -.098, 0, 0, .404, -.254, 0, .117, -.301, .662, .287, .414, 0, 0, .818, .243, 0, -4.450]';
coeffscontig2b=[.617,0,-.084,0,0,.369,-.554,0,.113,-.276,.756,.241,.427,0,0,0,.177,0,-3.509]';

%START THE SIMULATION
for n=1:10000
Countries=[Bahrain;Iran;Jordan;Kuwait;Qatar;Saudi;Syria;Turkey];
%Define a vector with binary inputs to record when states intervene
Interveners=zeros(8,1);
 
%START EACH RUN
for i=1:15
%Output predicted probabilities for each state into a vector
probs=normcdf(Countries*coeffscontig,0,1);
 
%Random draws to see who intervenes
draws=rand(8,1);
if draws(1)<probs(1)
    Interveners(1)=1;
end
    
if draws(2)<probs(2)
    Interveners(2)=1;
end
 
if draws(3)<probs(3)
    Interveners(3)=1;
end
 
if draws(4)<probs(4)
    Interveners(4)=1;
end
 
if draws(5)<probs(5)
    Interveners(5)=1;
end
 
if draws(6)<probs(6)
    Interveners(6)=1;
end
 
if draws(7)<probs(7)
    Interveners(7)=1;
end
 
if draws(8)<probs(8)
    Interveners(8)=1;
end
 
%Record progress
if i==5
   fiveyear(n)=sum(Interveners);
   Icounts(1,1)=Icounts(1,1)+Interveners(1);
   Icounts(2,1)=Icounts(2,1)+Interveners(2);
   Icounts(3,1)=Icounts(3,1)+Interveners(3);
   Icounts(4,1)=Icounts(4,1)+Interveners(4);
   Icounts(5,1)=Icounts(5,1)+Interveners(5);
   Icounts(6,1)=Icounts(6,1)+Interveners(6);
   Icounts(7,1)=Icounts(7,1)+Interveners(7);
   Icounts(8,1)=Icounts(8,1)+Interveners(8);
   if sum(Interveners)==0
       fivezero=fivezero+1;
   end
   if sum(Interveners)>0
       fiveone=fiveone+1;
   end
   if sum(Interveners)>1
       fivetwo=fivetwo+1;
   end
   if sum(Interveners)>2
       fivethree=fivethree+1;
   end
   if sum(Interveners)>3
       fivefour=fivefour+1;
   end   
end
 
if i==10
   tenyear(n)=sum(Interveners);
   Icounts(1,2)=Icounts(1,2)+Interveners(1);
   Icounts(2,2)=Icounts(2,2)+Interveners(2);
   Icounts(3,2)=Icounts(3,2)+Interveners(3);
   Icounts(4,2)=Icounts(4,2)+Interveners(4);
   Icounts(5,2)=Icounts(5,2)+Interveners(5);
   Icounts(6,2)=Icounts(6,2)+Interveners(6);
   Icounts(7,2)=Icounts(7,2)+Interveners(7);
   Icounts(8,2)=Icounts(8,2)+Interveners(8);
   if sum(Interveners)==0
       tenzero=tenzero+1;
   end
   if sum(Interveners)>0
       tenone=tenone+1;
   end
   if sum(Interveners)>1
       tentwo=tentwo+1;
   end
   if sum(Interveners)>2
       tenthree=tenthree+1;
   end
   if sum(Interveners)>3
       tenfour=tenfour+1;
   end
end
 
if i==15
   fifteenyear(n)=sum(Interveners);
   Icounts(1,3)=Icounts(1,3)+Interveners(1);
   Icounts(2,3)=Icounts(2,3)+Interveners(2);
   Icounts(3,3)=Icounts(3,3)+Interveners(3);
   Icounts(4,3)=Icounts(4,3)+Interveners(4);
   Icounts(5,3)=Icounts(5,3)+Interveners(5);
   Icounts(6,3)=Icounts(6,3)+Interveners(6);
   Icounts(7,3)=Icounts(7,3)+Interveners(7);
   Icounts(8,3)=Icounts(8,3)+Interveners(8);
   if sum(Interveners)==0
       fifzero=fifzero+1;
   end
   if sum(Interveners)>0
       fifone=fifone+1;
   end
   if sum(Interveners)>1
       fiftwo=fiftwo+1;
   end
   if sum(Interveners)>2
       fifthree=fifthree+1;
   end
   if sum(Interveners)>3
       fiffour=fiffour+1;
   end
end
 
%Advance arms race
Countries(1,6)=log(1.1*exp(Countries(1,6)*6.05))/6.05;
Countries(2,6)=log(1.1*exp(Countries(2,6)*6.05))/6.05;
Countries(3,6)=log(1.1*exp(Countries(3,6)*6.05))/6.05;
Countries(4,6)=log(1.1*exp(Countries(4,6)*6.05))/6.05;
Countries(5,6)=log(1.1*exp(Countries(5,6)*6.05))/6.05;
Countries(6,6)=log(1.1*exp(Countries(6,6)*6.05))/6.05;
Countries(7,6)=log(1.1*exp(Countries(7,6)*6.05))/6.05;
Countries(8,6)=log(1.1*exp(Countries(8,6)*6.05))/6.05;



end
end

%Show results
[Mean(fiveyear),std(fiveyear);Mean(tenyear),std(tenyear);Mean(fifteenyear),std(fifteenyear)]
Icounts/n
[fivezero,fiveone,fivetwo,fivethree,fivefour;tenzero,tenone,tentwo,tenthree,tenfour;fifzero,fifone,fiftwo,fifthree,fiffour]/n





%record prevint, multint
if sum(Interveners)>0
    Countries(1,17)=1;
    Countries(2,17)=1;
    Countries(3,17)=1;
    Countries(4,17)=1;
    Countries(5,17)=1;
    Countries(6,17)=1;
    Countries(7,17)=1;
    Countries(8,17)=1;
end
if sum(Interveners)>1
    Countries(1,18)=1;
    Countries(2,18)=1;
    Countries(3,18)=1;
    Countries(4,18)=1;
    Countries(5,18)=1;
    Countries(6,18)=1;
    Countries(7,18)=1;
    Countries(8,18)=1;
end
